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We improve the recently developed functional renormalization group (fRG) for 
impurities and boundaries in Luttinger liquids by including renormalization of the 
two-particle interaction, in addition to renormalization of the impurity potential. 
Explicit flow equations are derived for spinless lattice fermions with nearest neigh- 
bor interaction at zero temperature, and a fast algorithm for solving these equations 
for very large systems is presented. We compute spectral properties of single-particle 
excitations, and the oscillations in the density profile induced by impurities or bound- 
aries for chains with up to 10 6 lattice sites. The expected asymptotic power-laws at 
low energy or long distance are fully captured by the fRG. Results on the relevant 
energy scales and crossover phenomena at intermediate scales are also obtained. A 
comparison with numerical density matrix renormalization results for systems with 
up to 1000 sites shows that the fRG with the inclusion of vertex renormalization is 
remarkably accurate even for intermediate interaction strengths. 
PACS: 71.10.Pm, 72.10.-d, 73.21.Hb 



I. INTRODUCTION 

Low dimensional electron systems exhibit a rich variety of surprising effects which 
are due to the cooperative interplay of impurities and interactions. In one dimension 
even clean metallic systems are always strongly affected by interactions: at low energy 
scales physical properties obey anomalous power-laws, known as Luttinger liquid behav- 
ior, which is very different from the conventional Fermi liquid behavior describing most 
higher dimensional metals.^ 2 - In Luttinger liquids with repulsive interactions already a sin- 
gle static impurity is known to affect the low-energy properties drastically.^^™^ The 
impurity potential in a repulsive Luttinger liquid becomes dressed by long-range oscilla- 
tions which suppress the spectral weight for single-particle excitations near the impurity 
and also the conductance through the impurity down to zero in the low-energy limit. 

The asymptotic low-energy properties of Luttinger liquids with a single impurity are 
rather well understood. Universal power-laws and scaling functions have been obtained by 
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bosonization, conformal field theory, and exact solutions for the low-energy asymptotics 
in special integrable cases.-^ What remains to be developed, however, is a many-body 
method for microscopic models of interacting fermions with impurities, which does not 
only capture correctly the universal low-energy asymptotics, but allows one to compute 
observables on all energy scales, providing thus also non-universal properties, and in 
particular an answer to the important question at what scale the ultimate asymptotics 
sets in. That scale can indeed be surprisingly low, and the properties above it very 
different from the asymptotic behavior. Some of the non-universal properties can be 
computed numerically by the density matrix renormalization group (DMRG),— but this 
method is limited to lattice systems with about 1000 sites, and only a restricted set of 
observables can be evaluated with affordable computational effort. 

In the last few years it has been realized that the functional renormalization group 
(fRG) is a source of powerful new computation tools for interacting Fermi systems, es- 
pecially for low- dimensional systems with competing instabilities and entangled infrared 
singularities. The starting point of this approach is an exact hierarchy of differential flow 
equations for the Green or vertex functions of the system, which is obtained by taking 
derivatives with respect to an infrared cutoff A.— Approximations are then constructed 
by truncating the hierarchy and parametrizing the vertex functions with a manageable 
set of variables or functions.— 

A relatively simple fRG approximation for impurities in spinless Luttinger liquids 
has been developed recently by some of us.— The scheme starts from an fRG hierar- 
chy for one-particle irreducible vertex functions, as first derived in a field theoretical 
context by Wetterich 15 and Morris,— and for interacting Fermi systems by Salmhofer 
and Honerkamp. 17 A Matsubara frequency cutoff is used as the flow parameter. The hi- 
erarchy is then truncated already at first order, such that the 2-particle vertex remains 
unrenormalized, and the flow of the self-energy, which describes the renormalized impurity 
potential, is determined by the bare 2-particle interaction. No simplified parametrization 
of the self-energy is necessary such that the full spatial dependence of the renormalized 
impurity potential can be obtained for very large lattice systems. In spite of the striking 
simplicity of this scheme it was shown that the effects of a single static impurity in a 
spinless Luttinger liquid are fully captured qualitatively, and in the weak coupling limit 
also quantitatively— In particular, one obtains that impurity potentials in repulsive Lut- 
tinger liquids become effectively stronger at lower energy scales, and act ultimately as a 
weak link between two otherwise separate wires, as predicted by Kane and Fisher.^ The 
fRG also correctly yields the universal low energy power-laws with exponents that do 
not depend on the bare impurity strength. In addition, it was shown that the asymp- 
totic behavior holds typically only at very low energy scales and for very large systems, 
except for very strong bare impurities. 14 The fRG approach to impurities in Luttinger 
liquids was originally developed and tested for spectral densities of single-particle excita- 
tions, but has been applied very recently also to transport problems, such as persistent 
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currents in mesoscopic ring o 18 i 19 and the conductance of interacting wires connected to 
non-interacting leads . 19i2 Q The full power of the fRG with its ability to deal naturally also 
with complex crossover phenomena emerges most convincingly in the multi-scale problem 
posed by the transport through a resonant double barrier.— 

In the present work we further develop the fRG approach for impurities in Luttinger 
liquids by including 2-particle vertex renormalization, that is we go one step further in 
the hierarchy of flow equations than previously— For spinless fermions this extension 
does not matter qualitatively, as the lowest order is already qualitatively correct, but the 
quantitative accuracy of the results improves considerably in particular at intermediate 
interaction strengths. For spin-| systems, which we will treat in a subsequent work, 23 
vertex renormalization is necessary to take into account that backscattering of particles 
with opposite spins at opposite Fermi points scales to zero in the low energy limit. A 
crucial point is to devise an efficient parametrization of the vertex by a managable number 
of variables. Here we focus on spinless lattice fermions with nearest neighbor interaction 
as a prototype model. The bulk model is supplemented by site or hopping impurities. 
We also analyze the influence of boundaries, which can be viewed as infinite barriers or 
infinitesimal weak links. We choose to parametrize the vertex by a renormalized nearest 
neighbor interaction, which allows us to capture various advantageous features: the low 
energy flow of the vertex at kp in the pure system is obtained correctly to second order 
in the renormalized couplings; the non-universal contributions at finite energy scales are 
correct to second order in the bare interaction; the algorithm for the flow of the self-energy 
remains as fast as in the absence of vertex renormalization, such that one can easily deal 
with up to 10 6 lattice sites. We compute spectral properties of single-particle excitations 
near an impurity and the oscillations in the density profile induced by an impurity, or 
by a boundary. The accuracy of the calculation is checked by comparing with DMRG 
results for systems with up to 1000 lattice sites and with exact results for the asymptotic 
behavior, which can be obtained from the Bethe ansatz and bosonization. 

The article is structured as follows. In Sec. II we introduce the microscopic model 
and various types of impurities. The fRG formalism is developed in Sec. Ill, and worked 
out explicitly for the spinless fermion model with nearest neighbor interaction. Sec. IV is 
dedicated to results for the renormalized impurity potential, spectral properties, and the 
density profile. We finally conclude in Sec. V with an outline of promising extensions of 
the present work. 

II. MICROSCOPIC MODEL 

We consider the spinless fermion lattice model with nearest neighbor interaction and 
various types of impurity potentials. The Hamiltonian has the form 



H = H + Hj + H-, 



(1) 
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where the kinetic energy 

#o = -^(cj +1 c,+cjc J+1 ) (2) 

3 

is given by nearest neighbor hopping with an amplitude t, and 

Hi = U n 3 n 3+i ( 3 ) 



is a nearest neighbor interaction of strength U. We use standard second quantization 

"j and Cj 



notation, where cj and c, are creation and annihilation operators on site j, respectively, 



and rij = c] Cj is the local density operator. The impurity is represented by 

3,3' 

where Vj>j is a static potential. For "site impurities" 

V 3'3 = V 3 6 3f (•"») 

this potential is local. For the important special case of a single site impurity 

Vj = V6 Jjtt (6) 

the potential acts only on one site j . We also consider "hopping impurities" described 
by the non-local potential 

V 3'3 = V 3'f = -'j-.i ■ i $3' • ( 7 ) 
For the special case of a single hopping impurity, 

tj,j+i = (t " t) djj , (8) 

the hopping amplitude t is replaced by t' on the bond linking the sites j and j + 1. In 
the following we will set the bulk hopping amplitude t equal to one, that is all energies 
are given in units of t. 

The clean spinless fermion model Hq + Hi is exactly solvable via the Bethe ansatz.— 
The system is a Luttinger liquid for all particle densities n and any interaction strength 
except at half-filling for \U\ > 2. For U > 2 a charge density wave with wave vector 7r 
forms, for U < —2 the system undergoes phase separation. The Luttinger liquid parameter 
K p , which determines all the critical exponents of the liquid, can be computed exactly 
from the Bethe ansatz solution.— At half-filling K p is related to U by the simple explicit 
formula 

K; 1 = 1 arccos (-|) (9) 
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for \U\ < 2. 

We will compute physical properties for the above models on finite chains with open 
boundary conditions. The fRG equations can be solved easily for chains with up to 10 6 
sites, the DMRG can be applied for up to 10 3 sites. To avoid oscillations emerging from 
the ends of the chain we will sometimes attach non-interacting semi-infinite leads to the 
finite interacting chain, with a smoothly decaying interaction at the contact. 

III. FUNCTIONAL RENORMALIZATION GROUP 

This section is dedicated to a detailed presentation of the fRG equations for one- 
dimensional Fermi systems with impurities, including vertex renormalization, at zero 
temperature. We use the one-particle irreducible version of the functional RG . 1 ^ 1 ^ The 
starting point is an exact hierarchy of flow equations for the irreducible vertex functions, 
which is obtained by cutting off the infrared part of the free propagator on a scale A, 
and differentiating the generating functional for the vertex functions with respect to this 
scale. 1 '' 

A. Cutoff and flow equations 

The cutoff can be imposed in many different ways. The only requirement is that 
the infrared singularities must be regularized such that the flow equations allow for a 
regular perturbation expansion in powers of the renormalized two-particle vertex. For our 
purposes a sharp Matsubara frequency cutoff is the most efficient choice. A momentum 
cutoff is less suitable here since the impurity spoils momentum conservation. The cutoff is 
imposed by excluding modes with frequencies below scale A from the functional integral 
representation of the system, or equivalently, by introducing a regularized bare propagator 

G${iu) = 0(M - A)G (iu) . (10) 

Here Go is the bare propagator of the pure system, involving neither interactions nor 
impurities. Instead of the sharp cutoff imposed by the step function G one may also 
choose a smooth cutoff function, but the sharp cutoff has the advantage that it reduces 
the number of integration variables on the right hand side of the flow equations. Note 
that we will frequently write expressions which are well defined only if the sharp cutoff 
is viewed as a limit of increasingly sharp smooth cutoff functions. The suppression of 
frequencies below scale A affects all Green and vertex functions of the interacting system, 
which become thus functions of A. The original system is recovered in the limit A — ► 0. 

The first two equations in the exact hierarchy of flow equations for the one-particle 
irreducible m-particle vertex functions are depicted graphically in Figs. 1 and 2. The first 
equation yields the flow of the self-energy, which is related to the interacting propagator 
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by the usual Dyson equation 

G A = [(Gq) _1 — S A ]^ 1 . (11) 

Here and below G A , S A etc. are operators, which do not refer to any particular single- 
particle basis, unless we write matrix indices explicitly. The right hand side of the flow 
equation involves the 2-particle vertex T A , and the so-called single scale propagator S A 
defined as 

sA = GA N^r 1 ] GA = - rd^ d ~w rd^j • < 12 > 

Note that S A has support only on a single frequency scale: |a?| = A. The precise flow 
equation for the self-energy (sketched graphically in Fig. 1) reads 

As A (l', 1) = - - e^ 0+ S A (2, 2') r A (l', 2'; 1, 2) , (13) 

" 2,2' 

where (3 is the inverse temperature. The numbers 1, 2, etc. are a shorthand for Matsubara 
frequencies and labels for single-particle states such as site and spin indices. Note that 
u)i = u[ and u>2 = oj' 2 due to Matsubara frequency conservation. The exponential factor 
in the above equation is irrelevant at any finite A, but is necessary to define the initial 
conditions of the flow at A = A —>■ oo (see below). 

The right hand side of the flow equation for the 2-particle vertex T A , shown in Fig. 
2, involves T A itself, but also the 3-particle vertex T A . The latter could be computed 
from a flow equation involving the 4-particle vertex T A , and so on. To avoid this un- 
managable proliferation of vertex functions we make our first approximation: we neglect 
the contribution of the 3-particle vertex to the flow of T A . The coupled system of flow 
equations for T A and E A is then closed. The contribution of T A to T A is small as long 
as T A is sufficiently small, because T A is initially (at A ) zero and is generated only from 
terms of third order in T A .— A comparison of the RG results to exact DMRG results and 
exact scaling properties shows that the truncation error is often surprisingly small even 
for rather large interactions^ The explicit form of the truncated flow equation for the 
2-particle vertex reads 

^r A (l',2';l,2) = i^^G A (3,3')5 A (4,4') 

" 3,3' 4,4' 



X 



r A (l',2';3,4) r A (3 / ,4 / ;l,2) 
- r A (l', 4'; 1, 3) r A (3\ 2'; 4, 2) - (3 <-> 4, 3' <-> 4') 
+ T A (2', 4'; 1, 3) r A (3', 1'; 4, 2) + (3 <-> 4, 3' <-> 4') (14) 



The various contributions are shown diagrammatically in Fig. 3. 
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For a sharp frequency cutoff the frequency sums on the right hand side of the flow equa- 
tions can be carried out analytically in the zero temperature limit, where the Matsubara 
sums become integrals. At this point one has to deal with products of delta functions 
S(\uj\ — A) and expressions involving step functions 0(|u;| — A). These at first sight am- 
biguous expressions are well defined and unique if the sharp cutoff is implemented as a 
limit of increasingly sharp broadened cutoff functions O e , with the broadening parame- 
ter e tending to zero. The expressions can then be conveniently evaluated by using the 
following relation,— valid for arbitrary continuous functions /: 

6 e (x - A) f[Q e (x - A)] - 6{x - A) I f{t) dt , (15) 

Jo 

where 5 t = 0' e . Note that the functional form of e for finite e does not affect the result 
in the limit e — > 0. 

Instead of writing down the frequency-integrated flow equations in full generality, we 
first implement our second approximation: the frequency dependent flow of the renormal- 
ized 2-particle vertex T A is replaced by its value at vanishing (external) frequencies, such 
that T A remains frequency independent. As a consequence, also the self-energy becomes 
frequency independent. Since the bare interaction is frequency independent, neglecting 
the frequency dependence leads to errors only at second order (in the interaction strength) 
for the self-energy, and at third order for the vertex function at zero frequency. Besides 
the quantitative errors we miss qualitative properties related to the frequency dependence 
of the self-energy, such as the suppression of the one-particle spectral weight in the bulk 
of a pure Luttinger liquid. On the other hand, a comparison with exact numerical results 
and asymptotic analytical results shows that the impurity effects are not qualitatively 
affected by the frequency dependence of S, at least for weak interactions. 

Carrying out the frequency integration in the flow equation for the self-energy, and 
inserting a frequency independent 2-particle vertex, one obtains 

= E E e ^° + r^ ;1)2 , (16) 

w=±A 2,2' 

where the lower indices 1, 2, etc. label single-particle states (not frequencies) and 

G A (iuj) = [G^iiiu) -£ A ] _1 . (17) 

Note that G A has no jump at \u\ = A, in contrast to G A . The frequency integrated flow 
equation for the 2-particle vertex, evaluated at vanishing external frequencies, has the 
form 

A f a - 
dA L w* ~ 

2^ E EE [ 2 ^3,3' (^0 ^4,4' ( — r A ',2';3,4 ^3/ ,4/51,2 
w=±A 3,3' 4,4' 

+ ^3,3'(^) ^4,4' (^) ( — r A ',4';l,3 r A ',2';4,2 + l2',4';l,3 l3',l';4,2 
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The flow is determined uniquely by the differential flow equations and the initial condi- 
tions at A = oo. For A = oo the vertex functions are given by the bare interactions of the 
system. In particular, the flow of the 2-particle vertex starts from the antisymmetrized 
bare 2-particle interaction while m-particle vertices of higher order vanish at A = oo, in 
the absence of m-body interactions with m > 2. The self-energy at A = oo is given by 
the bare one-particle potential, that is by those one-particle terms which are not included 
already in Gq. We choose to include the translationally invariant bulk kinetic energy and 
the chemical potential in G , while the impurity (site or hopping) potential V is assigned 
to the initial condition of the self-energy. In a numerical solution the flow starts at some 
large finite initial cutoff Ao- Here one has to take into account that, due to the slow decay 
of the right hand side of the flow equation for E A at large A, the integration of the flow 
from A = oo to A = A yields a contribution which does not vanish in the limit A — > oo, 
but rather tends to a finite constant. Since G A 2 ,(iu;) ~^ $2,2'/ (^<-0 f° r M = A — > oo, this 
constant is easily determined as 

-- lim [ A ° dAYY e- 0+ I V 2 , x 2 = - Y I r 2 . x 2 , (19) 

J °° uj=±A 2,2' 2 

where Ii',2' ; i,2 is the bare antisymmetrized interaction. In summary, the initial conditions 
for the self-energy and the 2-particle vertex at A = A — > oo are 

S tl' = ^l,l' + ^A',2;l,2 (20) 
2 

ri°2';l,2 — h',2';l,2 , (21) 

where Vi t i> is the bare impurity potential. For the flow at A < Ao the factor e m0 in Eq. 
(16) can be discarded. 

B. Parametrization of 2-particle vertex 

For a finite lattice system with L sites the flow of the 2-particle vertex r A 2 ,. x 2 involves 
0(L 3 ) independent flowing variables, if translation invariance is assumed, and 0(L 4 ) 
variables, if the influence of the impurity on the flow of the 2-particle vertex is taken into 
account. For a treatment of large systems it is therefore necessary to reduce the number 
of variables by a suitable approximate parametrization of the vertex. In the low energy 
limit (small A) the flow is dominated by a very small number of variables anyway, the 
others being irrelevant according to standard RG arguments.— In particular, the frequency 
dependence of the vertex, discarded already above, is irrelevant for the flow of T A at small 
A. For larger A one can use perturbation theory as a guide for a simple but efficient 
parametrization of T A . 

We neglect the feedback of the self-energy into the flow of T A and replace G A by Go in 
Eq. (18). The renormalization of T A includes thus only bulk, not impurity, contributions, 
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such that T A remains translation invariant. While this turns out to be sufficient for 
capturing the effects of isolated impurities in otherwise pure systems, it is known that 
impurity contributions to vertex renormalization become important in macroscopically 
disordered systems.— For a more concrete treatment of the vertex renormalization, we 
have to focus on a specific model. In this article we restrict ourselves to the spinless 
fermion model with nearest neighbor interaction, leaving an extension to spin-| systems 
with more general interactions for subsequent work. 23 

For spinless fermions the 2-particle vertex and the self-energy are fully characterized 
by either site or momentum variables. In the low energy limit, the flow of the vertex is 
dominated by contributions with momenta close to the Fermi points, such that the right 
hand side of the flow equation is determined by momentum components of the vertex 
T k , k ,. ki k2 with k\, k2, k{, k' 2 = ±kp. Due to the antisymmetry of the vertex, there is only 
one such component which is non-zero: 

5 fA = ^k F ,~k F ;k F ,-k F ■ (22) 

In the low energy limit the momentum dependence of the vertex away from ±kp is ir- 
relevant. There are therefore many possible choices for the functional form of T A , k ,. ki k2 , 
which all lead to the correct low energy asymptotics. For a model with a bare nearest 
neighbor interaction U, a natural and efficient choice is to parametrize the flowing vertex 
simply by a renormalized nearest neighbor interaction U A , which leads to a real space 
vertex of the form 

r iU;^ = K,n - <Wi<W) ( 23 ) 

with C# = U A (6, 

ij'2-i + ^jij'2+1) • This yields the following structure in momentum 



space: 

n> A -MM = 2uA [ cos ( fc 'i - fc - cos (^ - h)} > ( 24 ) 

where the Kronecker delta implements momentum conservation (modulo 2tt). The flowing 
coupling constant U A is linked to the value of the vertex at the Fermi points by the relation 

g A = 2U A [1 - cos(2k F )} . (25) 
The flow equation for g A becomes 

^=s£/£< pp+pH+ra '> (26) 

uj=±A J 

with the particle-particle and particle-hole contributions 

PP = 2 ^p(^) G Q _ p {—iu) Tk F ,-k F ;p,-pTp-p;k F ,-k F 

PH = -[G° p (iuj)] 2 T kF ^. kF ^T A _ kF . p _ kF 

PH = Gp_ kp (iu) G p+kF (iuj) T A kp :P+ k F] k F ,p-k F ^p-k F ,k F ;p+k F ,-k F i (27) 
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where T A on the right hand side of the flow equation is given by Eq. (24). Using Eq. (25) 
to replace d\g A by d\U A on the left hand side of Eq. (26), one obtains a flow equation 
for U A of the simple form 

3 A U A = ft(A) (U A ) 2 . (28) 

The function ft(A) depends only on the cutoff A and the Fermi momentum kp. An explicit 
formula for h(A) can be obtained by carrying out the momentum integral in Eq. (26) using 
the residue theorem, as sketched in Appendix A. For finite systems the momentum integral 
should be replaced by a discrete momentum sum; however, this leads to sizable corrections 
only for very small systems. The flow equation (28) can be integrated to 

uA = T^wm • (29 > 

where H(A) is the primitive function of h(A) with H(A) —> for A — > oo. Integrating 
h(A) one obtains 



H(A) = -A + i Re 

Z7T 7T 



(4 - up A 2 - 2ifi (2 - /ig)A + /4 - 6/ig + 8 _ ijj^ _ t A - jfj, 
2(4-^)VA 2 -2^ A + 4-Atg 2 Sm 2 

(30) 



+ , ^tanh- 1 4 + ^ + ^ A 



2(4 -fJ$) 3 / 2 v /(4 + ^ + ZyUo A)2 + 4(A _2^ )2j 

where fio = —2coskp, while sinh -1 and tanh -1 denote the main branch of the inverse 
of the complex functions sinh and tanh, respectively. At half-filling, corresponding to 
kp = 7r/2, the function ft. (A) is particularly simple (see Appendix A), such that also U K 
has a relatively simple form: 

f/ A = -. . (31) 

1 + ( a - m) ^ 



in 



In Fig. 4 we show results for the renormalized nearest neighbor interaction U A as obtained 
from the flow equation at various densities n, for a bare interaction U — 1. While 
the renormalization does not follow any simple rule at intermediate scales A, all curves 
saturate at a finite value U* in the limit A — * 0, corresponding to a finite g*, as expected 
for a Luttinger liquid fixed point.— 

The above simplified flow equation yields not only the correct low energy asymptotics 
to second order in the renormalized vertex, but contains also all non-universal second 
order corrections to the vertex at ±kp at higher energy scales. Hence the resulting fixed 
point coupling g*, from which the Luttinger liquid parameter K p can be computed, is 
obtained correctly to second order in U. The calculation of K p , presented in Appendix 
B, illustrates how Luttinger liquid parameters can be related to the fixed point couplings 
obtained from the fRG. A comparison of the RG result for K p with exact results from the 
Bethe ansatz solution of the spinless fermion model shows that the vertex renormalization 
scheme described above is not only very simple, but also surprisingly accurate. 
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Parametrizing T A by a renormalized nearest neighbor interaction has the enormous 
advantage that the self-energy, as determined by the flow equation (16), is a tridiagonal 
matrix in real space, that is only the matrix elements E A ■ and E^ - ±1 are non-zero. Insert- 
ing r A from Eq. (23) into (16), one obtains the following simple coupled flow equations 
for the diagonal and off-diagonal matrix elements: 

dA^' j = ~ 27T Yl Gf +r , j+ r( iuj ) 

u=±Ar=±l 

^£<% ±1 (iu;). ( 32 ) 

Note that the self-energy enters also the right hand side of these equations, via G A = 
(Gq 1 — E A ) _1 . Since E A and Gq 1 are both tridiagonal in real space, the matrix inversion 
required to compute the diagonal and first off-diagonal elements of G A from E A can be 
carried out very efficiently, and thus for very large systems. In Appendix C we describe an 
algorithm for the numerical solution of the flow equation for E which scales only linearly 
with the size of the system. 



C. Observables 

In the next section we will present results for spectral properties of single-particle 
excitations near an impurity or boundary, and for the density profile. Here we describe 
how the relevant observables are computed from the solution of the flow equations. 



1. Single-particle excitations 

Integrating the flow equation for the self-energy E A down to A = yields the physical 
(cutoff-free) self-energy E, and the single-particle propagator G = (Gq 1 — E) -1 . From 
the latter spectral properties of single-particle excitations can be extracted. We focus on 
local spectral properties, which are described by the local spectral function 

pAu) = -- ImGjAuo + i0 + ) , (33) 

7T 

where Gjj(cu + i0 + ) is the local propagator, analytically continued to the real frequency 
axis from above. 

In our approximation the self-energy is frequency-independent and can therefore be 
viewed as an effective single-particle potential. The propagator G is thus the Green 
function of an effective single-particle Hamiltonian. In real space representation this 
Hamiltonian is given by the tridiagonal matrix h e s = ho + E, where the matrix elements 
of ho are the hopping amplitudes in H , Eq. (2). For a lattice with L sites this matrix has L 
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(including possible multiplicities) eigenvalues e\ and an orthonormal set of corresponding 
eigenvectors tp\. For the spectral function Pj(u>) one thus obtains a sum of ^-functions 

Pj (^) = Wx i S ( u ~&) ' (34) 

A 

where £,\ = t\ — fi, and the spectral weight w\j is the squared modulus of the amplitude 
of on site j. For large L the level spacing between neighboring eigenvalues is usually 
of order L -1 , except for one or a few levels outside the band edges which correspond to 
bound states. 

Due to even-odd effects etc. the spectral weight w\j generally varies quickly from 
one eigenvalue to the next one. A smooth function of u which suppresses these usually 
irrelevant finite-size details can be obtained by averaging over neighboring eigenvalues. 
In addition it is useful to divide by the level spacing between eigenvalues to obtain the 
conventional local density of states, which we denote by Dj(u). 



2. Density "profile. 

Boundaries or impurities induce a density profile with long-range Friedel oscillations, 
which are expected to decay with a power-law at long distances.-^ The expectation value 
of the local density rij could be computed from the local one-particle propagator Gjj, if 
G was known exactly. However, the approximate flow equations for E can be expected to 
describe the asymptotic behavior of G correctly only at long distances between creation 
and annihilation operator in time and/or space, while in the local density operator time 
and space variables coincide. In the standard RG terminology rij is a composite operator, 
which has to be renormalized separately.— 

To derive a flow equation for rij, we follow the usual procedure for the renormalization 
of correlation functions involving composite operators: we add a term (fijiij with a small 
field <j)j to the Hamiltonian, and take derivatives with respect to <j>j in the flow equations. 
The local density is given by 



(35) 

4> 3 =o 

where £l(4>j) is the grand canonical potential of the system in the presence of the field 4>j. 
Note that we use the same symbol rij for the density operator and its expectation value. 
In the presence of a cutoff A the grand canonical potential obeys the exact flow-equation 



OA (3 



l -Y. tr { e ^° + M" 1 ] [G A (iuj) - G^u)}} (36) 



in a short-hand matrix notation. This flow-equation follows from the functional flow 
equation for the generating functional for vertex functions^- and the relation between 
the grand canonical potential and the zero-particle vertex. At zero temperature the 
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Matsubara frequency sum becomes an integral which, for the sharp frequency cutoff (10), 
can be carried out analytically. This yields 



d_ 

OA 



QA = ^ E tr { e "° + ~ G ^ SA (^)]} • ( 37 ) 

oj=±A 



We attribute <pj to the "interaction" part of the Hamiltonian, not to H , such that Go 
remains independent of <pj. The self-energy is modified via the additional local and 
frequency-independent contribution <pj 5jj> to its initial value E^?, at scale Ao- 

The density profile can be obtained from the above equations and the flow equation 
for E A by computing the shift of Q generated by a small finite perturbation <pj, that is by 
numerical differentiation. Alternatively one may carry out the (^-derivative analytically 
in the flow equations, which yields a flow equation for the density in terms of the density 
response vertex. Taking the (fij derivative in (37) yields 

d 



^nf = -i- J2 ^ [e^ + G\iu) Rf(iu)] (38) 

uj=±A 

with the density response vertex 



d(j>j 



(39) 

4> 3 =o 



where the propagator G A is computed as previously, that is in the absence of 4>j- We 



compute the self-energy E in the presence of <f>j within the same approximation as pre- 
viously. It is thus determined from the flow equation (16) with a frequency-independent 
two-particle vertex T A . Taking a derivative of that equation with respect to <pj at <pj = 
yields the flow equation for the response vertex 



d 



4im = E EE^M43',3^(^)rf, 2 , ;1 , 2 . (40) 



OA 2tt 

w=±A 2,2' 3,3 



Note that Rf is frequency-independent in our approximation and that there is no contri- 



j 



bution from the ^-derivative of T , since we neglect self-energy contributions in the flow 
of the two-particle vertex. 

For spinless fermions with a (renormalized) nearest neighbor interaction, the matrix 
_R A is tridiagonal, that is only the components i? A ^ and Rf~n±i are non-zero, and their 
flow is given by 

qI R $;1,i = E EE E G trA iuJ ) R h',i'+r'G^ +r , tl+r (iuj) 

u;=±A I' r=±lr'=0,±l 

d U A 

gj^ R f;l,l±l = ^ E E E G t'( iu; ) R t;l',l'+r' G V+r',l±l( iu; ) ■ ( 41 ) 
ui=±A I' r'=0,±l 
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The right hand sides involve only one unrestricted lattice summation, and can be com- 
puted in O(L) time (see Appendix C). The initial condition for the response vertex is 
Rpui = SjiSw. The initial condition for the density is = |, for any filling, due to 
the slow convergence of the flow equation (38) at large frequencies, which yields a finite 
contribution to the integrated flow from A = oo to A for arbitrarily large finite A , as in 
the case of the self-energy discussed in more detail above. 

To avoid the interference of Friedel oscillations emerging from the impurity or one 
boundary with those coming from the (other) boundaries of our systems we suppress the 
influence of the latter by coupling the finite chain to semi-infinite non-interacting leads, 
with a smooth decay of the interaction at the contacts. The presence of leads modifies 
the self-energy in the interacting region only via a boundary term. In particular, a lead 
coupled to the first site of the interacting region via nearest neighbor hopping (with 
amplitude t — 1) yields the additional contribution 

to the first diagonal element of self-energy matrix in real space.— Note that this term 
is frequency dependent and independent of A. A lead coupled to the last site yields an 
analogous contribution to the last diagonal element of the self-energy. 



IV. RESULTS 



We now present and discuss explicit results for the self-energy, spectral properties 
near the impurity or boundary, and the density profile, as obtained from the fRG. A 
comparison with exact DMRG results is made for the spectral weight at the Fermi level 
and for the density profile. Some fRG results for the self-energy and spectral properties 
of the spinless fermion model have been published already previously, 1 " 4 but only for 
half-filling and without any vertex renormalization. Here we present results also away 
from half-filling and demonstrate the quantitative improvement of accuracy obtained by 
including vertex renormalization. In addition we present for the first time fRG results 
for the density profile. Using a faster algorithm (see Appendix C) than previously we are 
now able to solve the flow equations for much larger systems. 

The typical shape of the self-energy in the vicinity of an impurity can be seen in Fig. 5, 
where we plot the diagonal elements E 3 -j and the off-diagonal elements E/j+i near a site 
impurity of strength V = 1.5 added to the spinless fermion model with interaction strength 
U = 1 at quarter-filling. Recall that the self-energy is tridiagonal in real space and 
frequency independent within our treatment. The diagonal elements can be interpreted as 
a local effective potential, the off-diagonal elements as a non-local effective potential which 
renormalizes the hopping amplitudes. At long distances from the impurity both E^j and 
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Ejj+i tend to a constant. The former describes just a bulk shift of the chemical potential, 
the latter a bulk renormalization of the hopping amplitude toward larger values. The 
oscillations around the bulk shifts are generated by the impurity. The wave number of the 
oscillations is 2k p = n/2, where kp is the Fermi wave vector of the bulk system at quarter- 
filling. The amplitude of the oscillations of £ decays slower than the inverse distance 
from the impurity at intermediate length scales, but approaches a decay proportional to 
1/b ~~ Jo I f° r \j ~~ 3o\ ~^ 00 • This can be seen most clearly by plotting an effective exponent 
(3j for the decay of the oscillations, defined as the negative logarithmic derivative of the 
oscillation amplitude with respect to the distance \j — j \ . In Fig. 6 we show the effective 
exponent resulting from the oscillations of S^j as a function of the distance from a site 
impurity, for U = 1 and half-filling. The impurity is situated at the center of a long chain 
with L = 2 18 + 1 sites. To avoid interferences with oscillations from the boundaries we 
have attached semi-infinite non-interacting leads to the ends of the interacting chain, as 
described in the last paragraph of III.C.2. Only for relatively large impurity strengths 
the asymptotic regime corresponding to (3j = 1 is reached before finite size effects set 
in. For small V one can see that (5j increases from values below one, but the asymptotic 
long-distance behavior is cut off by the boundaries of the interacting region. For very 
small V (for example V = 0.01 in Fig. 6) we observe a plateau in f3j for intermediate 
distances from the impurity site. In this regime (3j is close to K p which can be understood 
by analytically solving our flow equations for small 1/.— 

The long-range 2/c^-oscillations of the self-energy lead to a marked suppression of the 
spectral weight for single-particle excitations at the Fermi level, that is at uj = 0. In Fig. 
7 we show the local density of states Dj(u) on the site next to a site impurity of strength 
V = 1.5 for the spinless fermion model at half-filling. The result for the interacting system 
at U = 1 is compared to the non-interacting case. Even-odd effects have been eliminated 
by averaging over neighboring eigenvalues (see Sec. III.C.l). 5-peaks outside the band 
edges corresponding to bound states are not plotted. The interaction leads to a global 
broadening of the band, which is due to an enhancement of the bulk hopping amplitude, 
and also to a strong suppression of Dj(uj) at low frequencies which is not present in the 
non-interacting system. For a finite system (here L = 1025) the spectral weight at the 
Fermi level remains finite, but tends to zero with increasing system size. In Fig. 8 we 
show results for the density of states choosing the same parameters as in Fig. 7, but now 
for densities away from half-filling: n = 1/4 and n = 3/4. In addition to the dip near 
uj = a second singularity appears at a finite frequency. This effect is due to the fact that 
a long-range potential with a wave number 2kp does not only strongly scatter states with 
momenta near kp, but also those with momenta close to 71 — kp. Indeed the singularity is 
situated at u — e^-kp — fi, where e& is the renormalized (bulk) dispersion. In the half- filled 
case only one singularity is seen simply because n — kp = kp for kp = tt/2. 

The spectral weight at the Fermi level is expected to vanish asymptotically as a power- 
law L~ aB with increasing system size, where as = K~ x — 1 is the boundary exponent 
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describing the power-law suppression of the density of states at the boundary of a semi- 
infinite chain. 7 That exponent depends only on the bulk parameters of the model, not 
on the impurity strength. For the spinless fermion model it can be computed exactly 
from the Bethe ansatz solution— We therefore analyze the asymptotic behavior of the 
spectral weight at the Fermi level by defining an effective exponent oi(L) as the negative 
logarithmic derivative of the spectral weight with respect to the system size, such that 
a(L) tends to a (positive) constant in case of a power- law suppression. In Fig. 9 we 
show results for a(L) as obtained from the fRG for the spinless fermion model at quarter- 
filling with up to about 10 6 sites, for a weak (U = 0.5) and an intermediate (U = 1.5) 
interaction parameter. The spectral weight has been computed either at a boundary, 
or near a hopping impurity of strength i! = 0.5. Results obtained from the RG without 
(upper panel) and with (lower panel) vertex renormalization are compared to exact DMRG 
results (for up to 512 sites) and the exact boundary exponents Ob, plotted as horizontal 
lines. The fRG results follow a power-law for large L, with the same asymptotic exponent 
for the boundary and impurity case, confirming thus the expected universality. However, 
the asymptotic regime is reached only for fairly large systems, even for the intermediate 
interaction strength U = 1.5. The comparison with the exact DMRG results and exact 
exponents shows that the fRG is also quantitatively rather accurate, and that the inclusion 
of vertex renormalization leads to a substantial improvement at intermediate coupling 
strength. Results for the effective exponent a in the case of a site impurity are shown in 
Fig. 10, at quarter-filling and for an interaction strength U — 1. The comparison of the 
different curves obtained for different impurity strengths confirms once again the expected 
asymptotic universality, and also how the asymptotic regime shifts rapidly toward larger 
systems as the bare impurity strength decreases. 

We finally discuss results for the density profile nj. Boundaries and impurities induce 
Friedel oscillations of the local density with a wave vector 2kp- In a non-interacting 
system these oscillations decay proportionally to the inverse distance from the boundary 
or impurity. In an interacting Luttinger liquid the Friedel oscillations are expected to 
decay as |j — jo I Kp at long distances |j — jo I- For a very weak impurity one expects a 
slower decay proportional to |j — jo\ 1 ~ 2Kp at intermediate distances, and a crossover to 
the asymptotic power-law with exponent K p at very long distances.— At intermediate 
distances the response of the density to a weak impurity can be treated in linear response 
theory, such that the density modulation is determined by the density- density response 
function at 2kp, which leads to the power-law decay with exponent 1K P — 1. In Fig. 11 we 
show fRG and DMRG results for the density profile rij for a spinless fermion chain with 
128 sites and interaction strength U — 1 at half-filling. The Friedel oscillations emerge 
from both boundaries and interfere in the center of the chain. The accuracy of the RG 
results is excellent for all j. For incommensurate filling factors the density profile looks 
more complicated. This can be seen in Fig. 12, where functional RG results are shown for 
the density modulation \rij — n\ near the boundary of a system with an average density 
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n = 0.393 and 8192 sites. For long distances from the boundary the oscillation amplitude 
has a well-defined envelope which fits to a power-law as a function of j. We now analyze 
the long-distance behavior of the amplitudes more closely for the half-filled case, and 
compare to exact results for the asymptotic exponents. In Fig. 13 we show fRG results 
for the amplitude of density oscillations emerging from an open boundary, for a very long 
spinless fermion chain with 2 19 + 1 sites and various interaction strengths U at half-filling. 
The other end of the chain (opposite to the open boundary) is smoothly connected to a 
non-interacting lead. In a log-log plot (upper panel of Fig. 13) the amplitude follows a 
straight line for almost all j, corresponding to a power-law dependence. Deviations from 
a perfect power-law can be seen more neatly by plotting the effective exponent aij, defined 
as the negative logarithmic derivative of the amplitude with respect to j (see the lower 
panel of Fig. 13). The effective exponent is almost constant except at very short distances 
or when j approaches the opposite end of the interacting chain, which is not surprising. 
From a comparison with the exact exponent (horizontal lines in the figure) one can assess 
the quantitative accuracy of the fRG results. Effective exponents describing the decay 
of Friedel oscillations generated by site impurities of various strengths are shown in Fig. 
14, for a half-filled spinless fermion chain with 2 18 + 1 sites and interaction U = 1. Both 
ends of the interacting chains are coupled to non-interacting leads to suppress oscillations 
coming from the boundaries. For strong impurities the results are close to the boundary 
result (cf. Fig. 13), as expected. For weaker impurities the oscillations decay more slowly, 
that is with a smaller exponent, and approach the boundary behavior only asymptotically 
at large distances (beyond the range of our chain for V < 1). For very weak impurities 
(V = 0.01 in Fig. 14) the oscillation amplitude follows a power-law corresponding to the 
linear response behavior with exponent 2K p — 1 at intermediate distances. We finally 
present some results for the effective exponent of the density oscillation decay in the case 
of an attractive interaction U = — 1, see Fig. 15. In that case the effective impurity 
strength should scale to zero at low energies and long distances.— Indeed, for weak and 
moderate bare impurity potentials the effective exponent in Fig. 15 approaches the linear 
density response exponent 2K p —l. Only for very strong impurities the density oscillations 
decay with the smaller boundary exponent over several length scales. 



V. CONCLUSION 

In summary, we have shown that the fRG provides an ideal tool for computing the 
intriguing properties of Luttinger liquids with static impurities. The method yields "ab 
initio" results for microscopic model systems at all energy scales from the Fermi energy 
to the ultimate low energy limit. We have demonstrated the power of the method by 
computing spectral properties of single-particle excitations as well as the oscillations in 
the density profile induced by impurities or boundaries for a spinless fermion model with 
nearest neighbor interaction. With the inclusion of vertex renormalization, in addition to 



18 



the renormalization of the effective impurity potential, our results agree remarkably well 
with exact asymptotic results and numerical DMRG data even for intermediate interaction 
strength. 

There is a broad range of interesting further applications and extensions of the fRG 
for impurities in Luttinger liquids beyond the scope of the present article: 

Spin-^ fermions: The inclusion of spin degrees of freedom requires a parametrization 
of the two-particle vertex with several coupling constants. We have already derived flow 
equations for the (extended) Hubbard model, with an effective vertex parametrized by 
local and nearest neighbor interactions.— 

Feedback of impurities on vertex: For isolated impurities the influence of impurities on 
the vertex renormalization is irrelevant for the asymptotic low energy or long distance 
behavior, although it may contribute quantitatively at intermediate scales. In the case of 
disordered systems with a finite impurity density the influence of the latter on the two- 
particle vertex is crucial and must be taken into account.^ In principle this is possible by 
computing the vertex flow with full propagators, which contain the renormalized impurity 
potential via the self-energy. 

Finite temperature: The fRG approach can be extended without major complications 
to finite temperature. This is particularly useful for studying the temperature dependence 
of transport properties. 30 

Transport: In the absence of inelastic processes (ImS = 0) the conductance of the 
interacting wire can be computed from the one-particle Green function in the presence of 
leads.- 3 ^ Several results have already been presented in short articles . 1 ?i 2 Vi 21 

Inelastic processes: Inelastic processes appear at second order in the interaction and 
can be included in the flow equations by inserting the second order vertex into the flow 
equation for the self-energy without neglecting its frequency dependence. This procedure 
would also capture the anomalous dimension of the bulk system, which is missing in the 
present work. 
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APPENDIX A: EVALUATION OF VERTEX FLOW 

Here we sketch some details concerning the explicit evaluation of the flow equations 
for the 2-particle vertex. Inserting the momentum structure of T A , (24), into the flow 
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equation (26) and replacing g by U on the left hand side yields 



auA - {uK? e r d £fM ( ai ) 



dA 27rshrfc F K J Q 2tt 



where 



2 sin 2 /cj? sin 2 p (cos kp — cosp) 2 [cos(2A;i?) — cosp] 2 

= & - W 11 ^ - z°p) ~ (™ - ^) 2 + - & kF )(™ - S3 ' 

Here £° = —2 cos k — fi with /x = —2 cos fci? is the bare dispersion relation relative to the 
bare Fermi level. Since f(p, u) can be written as a rational function of cosp and sinp, the 
p-integral can be carried out analytically using the substitution z = e ip and the residue 
theorem. After a lengthy but straightforward calculation one obtains the following result 
for the coefficient h(A) in (28): 



h(A) = -i- - Re 



i / 4 



2 yr ~" • Y * (/xo + ^A) 2 

3z/i^ - lOf/gA - 12i/i 2 (A 2 + 1) + 6A 3 /io + 18A// + 6iA 2 + iA 4 
X 7r(2/i + «A)(4 - /i 2 + A 2 - 2iA/x ) 2 

At half-filling (/i = 0), this reduces to 



(A3) 



A 



6 + A 



2 



(4 + A 2 ) 3 / 2 



(A4) 



APPENDIX B: CALCULATION OF K p 

In this appendix we show how the Luttinger liquid parameter K p , which determines the 
critical exponents of Luttinger liquids, can be computed from the fixed point couplings as 
obtained from the RG. A comparison of the RG result for K p with the exact Bethe ansatz 
result for the bulk model (without impurity) serves also as a check for the accuracy of 
our vertex renormalization. A relation between the fixed point couplings and K p can be 
established via the exact solution of the fixed point Hamiltonian of Luttinger liquids, the 
Luttinger model. 

For spinless fermions, K p is determined by the Luttinger model parameters g and Dp 

as 

= (bi) 

where g is the interaction between left and right movers and vf the effective Fermi velocity 
of the model, that is the slope of the (linear) dispersion relation, with a possible shift 
due to interactions between particles moving in the same direction (g 4 -coupling) already 
included.— We therefore need to extract g and vp from the RG flow in the limit A — > 0. 
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In order to obtain K p correctly to order U 2 , it is sufficient to obtain vp correctly to linear 
order in U. 

The Luttinger model interaction g and the fixed point coupling g* = r£l*° fc fc _ kp 
from the RG are not simply identical, in contrast to what one might naively expect. To 
find the true relation between g and g*, one has to take into account that the forward 
scattering limit of the dynamical 2-particle vertex is generally not unique (in the absence 
of cutoffs), and depends on whether momentum or frequency transfers tend to zero first. 
This ambiguity is well-known in Fermi liquid theory, where it leads to the distinction 
between quasi-particle interactions and scattering amplitudes, 27 but is equally present in 
Luttinger liquids, for the same reason in all cases: the ambiguity of the small momentum, 
small frequency limit of particle-hole propagators contributing to the vertex function. In 
the dynamical limit, where the momentum transfer q vanishes first, the singular particle- 
hole propagators do not contribute. In Fermi liquids this limit yields the quasi-particle 
interaction. In the opposite static limit the frequency transfer v vanishes first and particle- 
hole propagators yield a finite contribution. In the presence of an infrared cutoff A > 
the forward scattering limit of the vertex function is unique, since the ambiguity in the 
particle-hole propagator is due to the infrared pole of the single-particle propagator. Hence 
V k F ,-k F -,k F ,-k F is wel1 defined. However, V^ F _ kF . kF _ kF and also its limit for A -> depend 
on the choice for the cutoff function. For a momentum cutoff, which excludes states with 
excitation energies below A around the Fermi points, particle-hole excitations with small 
momentum transfers q are impossible. Hence particle-hole propagators with infinitesimal 
q do not contribute to the vertex at any A > 0, such that T kp _ kp . kF _ kp converges to the 
dynamical forward scattering limit, which is simply given by the bare coupling constant g 
in the Luttinger model. For a frequency cutoff the particle-hole propagators with vanishing 
momentum and frequency transfer yield a finite contribution at A > 0, which tends to 
the static limit for A — > 0. This can be seen directly by integrating XL=±A I [Gp(iuj)] 2 
over A from infinity to zero. Hence the vertex T kp _ kp . kF _ kp obtained from our frequency 
cutoff RG tends to the static forward scattering limit. 

For the Luttinger model, the static forward scattering limit of the vertex can be ob- 
tained from the dynamical effective interaction between left and right movers D(q,iv), 
which is defined as the sum of particle-hole chains 

D(q, iv)=g + g n°_ (q, iv) gU° + (q,zu) g + ■ ■ ■ = - 9 . , (B2) 

1 — g z lL_[q, iv) yl+{q, w) 

where 

Ul(q,iv) = ±- q - (B3) 

is the bare particle-hole bubble for right (+) and left (-) movers. Note that only odd 
powers of g contribute to the effective interaction between left and right movers. This 
effective interaction appears naturally in the exact solution of the Luttinger model via 
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— TTVp 



+ V(™f) 2 + (g*) 2 } (B5) 



Ward identities.— For the static limit one obtains 

y<«-°> = iHsfai. (B4) 

which we identify with our fixed point coupling g* as obtained from the RG with frequency 
cutoff. Inverting this relation between g and g* we obtain 

27TVF 

9* 

For spinless fermions the difference between g and g* appears only at third order in the 
coupling, but for models with spin the distinction becomes important already at second 
order. 

The Fermi velocity vf can be computed from the (frequency independent) self-energy 
in momentum space as 

v F = v° F +dj: k \ kF , (B6) 

where v F = dk€k\k F is the bare Fermi velocity. The self-energy is computed from the flow 
equation (32), which can be rewritten in momentum space as 

d A U A ^ f dp 1 - cos(/c - p) 



where £ p = e p — //. The chemical potential /x has to be fixed by the final condition 
£k F + Efc F = 0, where fcp = vrn depends only on the density, not the interaction. From 
the tridiagonal structure of E in real space, but also from the above expression it follows 
that E A has the form E A = a A + b A cos A;. The functional flow equation for yields a 
coupled set of ordinary flow equations for the coefficients a A and b A , with initial conditions 
a A ° = U and b A ° = 0. The momentum integrals can be evaluated analytically via the 
residue theorem, such that the remaining set of two coupled differential equations (with 
U A as input) can be easily solved numerically. The result for vp is correct at least to 
first order in U, but not necessarily to second order, since our simplified parametrization 
captures the 2-particle vertex correctly to second order only at the Fermi points. 

Inserting g and vf into the Luttinger model formula (Bl) we can now compute K p as a 
function of U and density for the microscopic spinless fermion model. In Fig. 16 we show 
results for K p (U) for various fixed densities, which are compared to the exact result from 
the Bethe ansatz solution^ in the inset. The RG results are indeed correct to order U 2 
for small U, as expected, and they are surprisingly accurate also for larger values of U. 



APPENDIX C: NUMERICAL SOLUTION OF S FLOW 



For A < Aq < oo the flow equation for the self-energy (16) can be written as 



— e a -_iVr A 



2 Re 



(CI) 



2,2' 
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In order to compute its right-hand side, one needs to invert the tridiagonal matrix 



T = G 1 (zA)-£ / 



(C2) 



where T is complex symmetric (not hermitean) with diagonal elements aj := zA + /i — E^, 
i — 1, . . . , L, and first off-diagonal elements bi :— t — i — 1, . . . , L — 1. Note that 

Im(aj) = A > such that T is non-singular and its inverse well-defined. 

The inverse G A (iA) = T -1 is not tridiagonal but a full matrix which can be computed 
by standard methods in 0(L 2 ) time. However, for an interaction that does not extend 
beyond nearest neighbors on the lattice, only the tridiagonal part of G is required, which 
can be computed in O(L) time, such that much larger lattices can be treated. We shall 
first explain how this is done and then present the resulting algorithm that can directly 
be incorporated into a computer program. 

Under certain assumptions (see below), a matrix can be uniquely factorized into a 
lower unit triangular matrix L, a diagonal matrix D, and an upper unit triangular matrix 
U ("LDU factorisation"): T = L ■ D ■ U.- 2 For a tridiagonal matrix T the unit triangular 
matrices L and U are in fact unit bidiagonal: their matrix elements are one on the 
diagonal, and only the first off-diagonal is nonzero. Since our T is symmetric we have 
L = U T . Thus we obtain a factorization of the form 



T 



U +T D + U- 
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\ 
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D2 



\ f 1 ut 
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D: 
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7 



where the label "+" distinguishes this factorization from another one used below. The 
prescription to compute the elements Df and Uf is well known and can be found for 
example in Ref.^. Starting in the upper left corner one proceeds to increasing row and 
column numbers until one arrives at the lower right corner of T: 

D+:=a 1} U?:=bi/D+, D+ +1 := a l+1 - b t U+ (i = 1, . . . , L - 1) . (C3) 

This works well since in our case Im(D^) > A > 0, such that one never divides by zero. 

To compute the inverse G = T , one could directly calculate (U + y 1 (D + )- 1 (U +T )~ 1 . 
It is however easier and more accurate to find the inverse by solving the linear system of 
equations TG = 1, where 1 is the identity matrix, by "backsubstitution" . To be specific, 
consider the i th column vector G. , of G: 



Gi = TG. ti = U +T (D+U+G.,) = U +T 9l , U + G., = (D 



(C4) 



where is the i unit vector. The first step is to solve the linear system U +T = ei 
for gi, and the second step to solve U + G. t i = (D + )~ 1 g i for G./ t . To solve a tridiagonal 
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linear system for one vector takes O(L) time, so solving for the full inverse matrix G takes 
0(L 2 ) time. 

Now we shall derive an algorithm to compute the elements of gi and G.^. Begin with 
the last column i = L: U +T gL = &l can be solved from the first to the last row and 
gives gi, = e^. Next U + G.^ = (D + )~ 1 ei J can be solved starting from the last row, 
Gt,L — V^l- F rom there one can work upwards by backsubstitution, Gj^ = —Uj~Gj+i t L 
(j = 1, . . . , L — 1). For the other columns i < L, one cannot take the shortcut and has 
to solve both linear systems for and G.^. But it is now important to realize that for 
any column vector G.^+i, if we somehow know the diagonal element Gj+i^+i, the next 
element above the diagonal is 

G i>i+1 = -U?G i+1 s +1 (i = l,...,L-l). (C5) 

Thus we have a prescription how to go up one row in G. Together with the symmetry 
of G, that is G^i+i = Gj+i^, which follows from the symmetry of T, we get the first off- 
diagonal element one column to the left without solving the two linear systems in (jC4|) . 
So it is possible to compute directly the tridiagonal part of the inverse. However, there is 
another algorithm which is much more accurate for near-singular matrices at the end of 
the RG flow: the double factorization.— It does not rely on the symmetry of G but uses 
the complementary "UDL" factorization 

T = U~D~L- = U-D-U- T , (C6) 

where the matrix elements are obtained as 

D L :=a L , Ur-^bi/D^, D7 : = a, - b { U~ (i = L - 1, . . . , 1) . (C7) 

We proceed as for the LDU factorization above and get 

Gi t i = 1/L>r (C8) 
G hl+1 = -Ur ■ G hl . (C9) 

We can combine equations (jC5|) and (jC9|) to relate consecutive diagonal elements: 

G i+W = -G hl+1 /Ut = G iti ■ Ur/U+ = G iti ■ Dt/Dr +l . (CIO) 

Thus we start with (|C8|) and use the U~D~L~ decomposition to go one matrix element 
to the right in the inverse matrix, from the diagonal to the first off-diagonal (|C9J) . while 
the L + D + U + decomposition allows to go down by one, back to the next diagonal element 
(|C5|) . There is no need to compute the full inverse matrix. 

One can implement the algorithm without knowing the derivation by using equations 
(jC3j) and (jC7j) - (|C10|) . One can further eliminate the U's using equations (jC3j) and 
()C7|) and implement the algorithm such that only the input vectors a i: bi and the output 
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vectors G^, (5^+1 enter the temporary storage. This double factorization is numerically 
accurate to more than 10 significant digits (using double precision) even for large lattices 
(10 6 sites) and almost singular matrices with |ai| ~ 10~ 15 which appear at the end of the 
flow for half filling. 

The right hand side of the flow equations (41) for the density response vertex R A can 
be computed in O(L) time using the fact that the upper triangular part of the inverse 
of a tridiagonal matrix (here G) is the upper triangular part of the outer product of two 
vectors.—** 
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FIG. 2: Flow equation for the 2-particle vertex T A 
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FIG. 3: Flow equation for the 2-particle vertex at 1-loop level with particle-particle and particle- 
hole channels written explicitly. 
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FIG. 4: Flow of the renormalized nearest neighbor interaction C/ A for the spinless fermion model, 
for U = 1 and various densities n. 
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FIG. 5: Self-energy near a site impurity of strength V = 1.5 for the spinless fermion model at 
quarter-filling and interaction strength U = 1; the impurity is situated at the center of a chain 
with L = 1025 sites. 
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FIG. 6: Effective exponent for the decay of oscillations of Y,jj as a function of the distance 
from a site impurity of strengths V = 0.01, 0.1, 0.3, 1, 10 (from bottom to top), for the spinless 
fermion model at half-filling and interaction strength £7 = 1; the impurity is situated at the 
center of a chain with L = 2 18 + 1 sites. 
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FIG. 7: Local density of states on the site next to a site impurity of strength V = 1.5 for spinless 
fermions at half-filling and U = 1; the impurity is situated at the center of a chain with L = 1025 
sites; the non-interacting case U = is shown for comparison. 




FIG. 8: Local density of states on the site next to a site impurity as in Fig. 7 (same parameters), 
but now for densities n = 1/4 and n = 3/4. 
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FIG. 9: Logarithmic derivative of the spectral weight at the Fermi level near a boundary (solid 
lines) or hopping impurity (dashed lines) as a function of system size L, for spinless fermions 
at quarter-filling and interaction strength U = 0.5 (circles) or U = 1.5 (squares); upper panel: 
without vertex renormalization, lower panel: with vertex renormalization; the open symbols 
are fRG, the filled symbols DMRG results; the horizontal lines represent the exact boundary 
exponents for U = 0.5 and U = 1.5. In the boundary case (solid lines) the spectral weight has 
been taken on the first site of a homogeneous chain, in the impurity case (dashed lines) on one 
of the two sites next to a hopping impurity t' = 0.5 in the center of the chain. 
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FIG. 10: Logarithmic derivative of the spectral weight at the Fermi level near a boundary (solid 
line) or near a site impurity (dashed lines) as a function of system size L, for spinless fermions 
at quarter-filling and interaction strength U = 1. In the boundary case the spectral weight has 
been taken on the first site of a homogeneous chain, in the impurity case on the site next to a 
site impurity of strength V in the center of the chain. The horizontal line represents the exact 
boundary exponent for U = 1. 




FIG. 11: Density profile rij for a spinless fermion chain with 128 sites and interaction strength 
U = 1 at half-filling; fRG results are compared to exact DMRG results. 
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FIG. 12: Density modulation \rij — n\ as a function of the distance from a boundary, for spinless 
fermions with interaction strength U = 1 and (average) density n = 0.393 on a chain with 8192 
sites; the broken straight line is a power-law fit to the envelope of the oscillation amplitudes. 



33 



' — - 



io- 1 
io- 2 

10~ 3 

io~ 4 

1(T 6 

io- 7 



"A. 

U3. 



10" 



*^A X>-^ 



IO 1 



s. 



B-. 



ex. 



Ek 



•o' 



"~ex 



--O-- 
--A-- 
-43-- 

I L 



U = 0.1 

U = 0.5 

U = 1.0 

U = 1.5 
i . j_ 



10 2 



10 3 



10 4 



10 5 



J 



s 

s 
o 



1 

0.9 
0.8 



1 1 — n 1 1 — 1-| 1 1 — 1-| 1 1 — 1-| 1— 

©- -e-- e— e— e— ©-©-- e— ©-- e— e— e --e---©- -e- -e- -e- -° 

Al 



A - A - A- -A"- -A- A"-" A-"-" A- &--A--&-S - A -A - A^ -A- 




FIG. 13: Amplitude (envelope) of oscillations of the density profile rij induced by a boundary 
as a function of the distance from the boundary, for spinless fermions with various interaction 
strengths U at half-filling; the interacting chain with 2 19 + 1 sites is coupled to a semi-infinite 
non-interacting lead at one end (opposite to the boundary); upper panel: log-log plot of the 
amplitude, lower panel: effective exponents for the decay, and the exact asymptotic exponents 
as horizontal lines. 
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FIG. 14: Effective exponent for the decay of density oscillations as a function of the distance 
from a site impurity of strengths V = 0.01, 0.1, 0.3, 1, 3, 10 (from bottom to top); the impurity 
is situated at the center of a spinless fermion chain with 2 18 + 1 sites and interaction strength 
U = 1 at half-filling; the interacting chain is coupled to semi-infinite non-interacting leads at 
both ends. 
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FIG. 15: Effective exponent for the decay of density oscillations as a function of the distance 
from a site impurity of strengths V = 0.1, 1, 10, 100, 1000 (from top to bottom) for the same 
chain as in Fig. 14 but now with an attractive interaction, U = — 1. 
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FIG. 16: Luttinger liquid parameter K p as a function of U at various densities (as in Fig. 4) for 
the spinless fermion model; the inset shows the difference between the RG result and the exact 
Bethe ansatz result for K p . 



